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1.  Introduction 

Composite  materials  are  widely  used  in  aerospace  and  automo¬ 
bile  structures,  sporting  goods  and  military  equipment  since  they 
show  superior  specific  stiffness  and  strength.  Delamination  is  one 
of  the  most  significant  damage  modes  in  laminated  composites 
and  is  a  key  element  of  consideration  in  design  and  operation.  It 
is  well  known  that  translaminar  reinforcements  can  enhance  the 
delamination  resistance  and  damage  tolerance  of  composites  sig¬ 
nificantly  [1-4].  It  is  also  known  that  sparse  through-the-thickness 
reinforcements  are  not  beneficial,  and  at  the  same  time  excessive 
amount  of  reinforcement  is  also  deleterious  as  it  results  in  damage 
to  the  composite  structure  and  degradation  in  its  properties  [4].  The 
effect  of  z-fiber  on  delamination  of  composites  has  been  extensively 
studied.  Experimental  approaches  involve  measurement  of  appar¬ 
ent  fracture  toughness  [5],  characterization  of  pull-out  process  of 
z-fiber  [6,7]  and  crack  propagation  [7].  Analytic  models  are  also 
available  to  find  apparent  fracture  toughness  and  to  establish  the 
relationship  among  various  design  parameters  [8-13].  Most  ana¬ 
lytic  models  used  beam  or  plate  models  in  which  stitches  are  repre¬ 
sented  by  an  effective  bridging  force.  Either  Euler-Bernoulli  beams 
[9-13]  or  Timoshenko  beams  [8,9]  are  used  in  association  with 
discrete  element  [10]  or  continuous  element  [8,9,11-13]  for  the 
stitches.  Numerical  studies  focused  on  crack  initiation  and  progres¬ 
sive  delamination  in  which  J-integral  [7],  virtual  crack  closure 
technique  (VCCT)  [15]  and  cohesive  zone  method  [14,16]  were 
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primarily  used  along  with  beam  elements  [14],  plane  strain 
elements  [7],  solid  elements  [16]  or  shell  elements  [15].  Sankar 
and  Hu  [17]  used  beam  elements  and  spring  elements  similar  to 
cohesive  zone  elements  to  simulate  dynamic  delamination  propa¬ 
gation  due  to  impact. 

While  FE  models  that  account  for  each  and  every  reinforcement 
may  be  more  accurate  in  predicting  the  damage  characteristics, 
they  tend  to  be  computationally  expensive  for  realistic  structures. 
On  the  other  hand  analytical  models  that  smear  the  reinforcements 
as  a  continuous  element  are  efficient  and  provide  insight  into  the 
mechanics,  but  tend  to  be  unrealistic  in  some  situations. 

A  summary  of  past  work  in  this  area  is  provided  in  Table  1 .  In  the 
following  we  briefly  describe  the  salient  features  of  the  work  listed 
in  that  table.  Cartie  [7]  suggested  a  bi-linear  bridging  law  through 
experiments  and  conducted  FE  simulation  using  plain  strain  ele¬ 
ments.  Robinson  and  Das  [11],  Mabson  and  Deobald  [12]  and  Byrd 
and  Birman  [13]  used  linear  softening  bridging  model  whereas 
Ratcliffe  and  O’brien  [10]  proposed  discrete  element  analytical 
model.  For  delamination  modeling  of  z-pinned  DCB,  Dantuluri  et 
al.  [14]  used  2D  cohesive  elements  incorporating  beam  elements 
and  a  bi-linear  bridging  law.  Ratcliffe  and  Krueger  [16]  used  3D 
cohesive  model  with  solid  elements.  Grassi  and  Zhang  [15]  used 
VCCT  for  delamination  modeling.  It  is  clear  from  the  table  that  var¬ 
ious  researchers  have  used  different  approaches  in  modeling  the 
stated  problem.  Our  goal  here  is  to  use  both  analytical  and  numer¬ 
ical  approaches  and  compare  them  with  available  experimental 
data.  In  particular,  we  have  developed  a  non-dimensional  model 
that  will  be  useful  in  the  design  of  translaminar  reinforcements 
for  a  given  application.  In  particular  we  have  used  a  3D  shell  model 
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Table  1 

Summary  of  results  obtained  by  various  researchers  and  present  study. 


Ref. 

z-Pin  in  analytic  model 

z-Pin  model  in  FEM 

Crack  propagation 
criterion  in  FEM 

Sub-laminate 
model  in  FEM 

Comparison  with 
experimental  data 

Cartie  [7] 

_ 

Bi-linear  traction 

J  integral 

Plane  strain 

Yes 

Mabson  and  Deobald 
[12] 

Distributed  force 
(non-dimensional  form-1) 

- 

- 

- 

- 

Byrd  and  Birman  [13] 

Distributed  force 
(dimensional  form-1) 

- 

- 

- 

- 

Robinson  and  Das  [11] 

Distributed  force 
(dimensional  form-2) 

- 

Yes 

Grassi  and  Zhang  [15] 

- 

Nonlinear  spring  (bi-linear  function) 

VCCT 

Shell  element 

Yes 

Ratcliffe  and  O’brien 

Discrete  force 

- 

- 

- 

Yes 

[10] 

Dantuluri 

- 

Nonlinear  spring  (bi-linear  function) 

2D  cohesive  element 

Beam  element 

Yes 

et  al.  [14] 

Ratcliffe  and  Krueger 

- 

Nonlinear  spring  (bi-linear  function) 

3D  cohesive  element 

Solid  element 

- 

[16] 

Present  study 

Distributed  force 
(non-dimensional  form-2) 

Nonlinear  spring  (linear  softening 
function  and  tri-linear  bridging  law) 

3D  cohesive  element 

Shell  element 

Yes 

in  conjunction  with  3D  cohesive  elements.  Both  linear  softening 
law  and  more  realistic  tri-linear  bridging  law  are  used  to  model 
the  z-pins.  In  order  to  clarify  the  role  of  z-pins,  a  non-dimensional 
analytical  model  is  proposed.  The  solution  of  the  non-dimensional 
equation  uses  an  iterative  procedure.  We  have  also  derived  an 
expression  for  the  maximum  density  of  z-pins  that  can  be  allowed 
before  the  beam  fails  otherwise.  The  efficacy  of  the  analytical  model 
is  verified  by  finite  element  simulation  of  the  DCB  specimen.  We 
have  used  the  example  given  in  [7]  for  this  purpose.  In  the  FE  sim¬ 
ulation,  the  ligaments  of  the  DCB  are  modeled  using  shell  elements. 
Cohesive  elements  are  used  to  simulate  the  delamination  and  dis¬ 
crete  nonlinear  elements  are  used  to  model  the  z-pins.  The  agree¬ 
ment  between  the  analytical  model  and  FE  simulations  is  found 
to  be  excellent  for  various  results  such  as  load-deflection,  load- 
crack  length,  and  effective  fracture  toughness. 

It  was  found  that  steady-state  bridging  length  and  maximum 
apparent  fracture  toughness  can  be  related  to  interlaminar  fracture 
toughness  of  the  composite  and  maximum  frictional  force  sup¬ 
ported  by  z-pins.  The  relationships  among  various  parameters  are 
non-dimensionalized  and  the  maximum  pin  friction  that  can  be 


allowed  before  the  composite  beam  itself  fails  is  calculated.  The 
non-dimensional  analytical  model  could  be  a  useful  design  tool 
in  selecting  z-pins  for  composite  structures  to  improve  interlami¬ 
nar  fracture  toughness. 

2.  Non-dimensional  analytical  model 

Consider  a  z-pinned  composite  DCB  specimen  of  thickness  2 h 
with  initial  crack  length  of  a0  and  subjected  to  Mode  I  loading  as 
shown  in  Fig.  la.  A  pair  of  transverse  forces  F  is  applied  at  the  tip 
of  the  DCB.  Due  to  the  applied  transverse  force  F  the  initial  crack 
tends  to  reach  the  current  crack  length  denoted  by  a  (Fig.  lb).  When 
the  crack  reaches  the  region  reinforced  by  z-pins,  a  bridging  zone  of 
length  c  begins  to  develop  (Fig.  lc).  We  can  define  an  apparent  crack 
length  ap  which  is  the  length  of  the  crack  up  to  the  beginning  of 
bridging  zone.  Current  crack  length  (a)  is  the  sum  of  the  apparent 
crack  length  (ap)  and  the  bridging  length  (c).  When  the  z-pins  start 
to  be  completely  pulled  out  of  the  composite  as  shown  in  Fig.  Id, 
the  bridging  zone  becomes  fully  developed  and  a  new  fracture  sur¬ 
face  is  created  in  the  z-pinned  zone.  Therefore  in  the  bridging  zone 
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Fig.  1.  Schematic  of  initial  geometry  and  development  of  bridging  zone  by  z-pins  in  the  composite  DCB. 
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the  pull-out  of  pins  is  partial.  Beyond  the  current  crack  tip  the  pins 
are  assumed  to  be  intact.  The  relative  deflection  at  the  tip  of  the  DCB 
is  denoted  by  3  as  shown  in  Fig.  le.  Our  goal  is  to  determine  the 
relationships  among  F,  a,  c,  and  3.  Once  the  bridging  zone  length  c 
is  determined  other  parameters  such  as  F,  a,  and  3  can  be  found 
by  solving  the  governing  equations.  In  the  following  we  describe 
the  non-dimensional  equation  governing  the  bridging  zone  and 
procedures  to  determine  the  bridging  length  for  a  given  load. 

We  assume  that  the  shear  deformation  is  negligible  and  use  Eu- 
ler-Bernoulli  beam  equations  to  model  the  ligaments  of  the  DCB  as 
shown  in  Fig.  2b.  We  also  assume  that  z-pins  are  rigid  and  there  is 
sufficient  friction  between  the  z-pins  and  the  composite  material 
surrounding  it.  For  simplicity,  the  relation  between  the  frictional 
force  (f)  and  pull  out  or  slip  distance  (ds)  is  idealized  as  a  linear  soft¬ 
ening  function  as  shown  in  Fig.  3.  The  validity  of  this  assumption 
will  be  later  verified  in  FE  simulation.  When  the  pins  are  intact  they 
can  exert  a  maximum  friction  force  of/m.  As  the  pins  pull  out  of  the 
material,  the  loss  of  friction  is  proportional  to  the  pullout  distance 
ds.  When  the  pin  is  completely  pulled  out  of  the  beam,  the  friction 
force  reduces  to  zero.  Thus  th ef-ds  relationship  is  given  by: 

/=/m(1-|),  (1) 

where  h  is  half  the  thickness  of  DCB. 

Although  the  resistance  offered  by  the  pins  on  the  beam  is  dis¬ 
creet,  for  the  purpose  of  the  analytical  model  we  smear  the  discrete 
pin  resistance  as  continuous  distributed  traction  (p)  acting  on  the 
crack  surfaces  as  shown  in  Fig.  2c.  Then  the  traction  can  be  derived 
as: 


(a) 


Fig.  2.  Idealization  of  bridging  zone  using  a  beam  model. 


P  =  Nf  =  Nfm(  1 


(2) 


where  N  is  z-pin  density  expressed  as  number  of  z-pins  per  unit 
area. 

The  Euler-Bernoulli  beam  equation  for  one  of  the  ligaments,  say 
upper  ligament  (Fig.  2c),  of  the  DCB  can  be  written  as: 


E,dw 

El-rr  =  -bp 
dx 


(3) 


where  b  is  the  beam  width.  The  effective  bending  rigidity  of  one  of 
the  ligaments  of  the  DCB  is  represented  by  the  term  EL  Note  that  the 
flexural  rigidity  for  a  laminated  composite  can  be  taken  from  the 
bending  stiffness  matrix  in  classical  plate  laminate  theory.  Substi¬ 
tuting  for  p  from  Eq.  (2)  into  Eq.  (3),  we  obtain: 


_,d4w  ds 

E  ~dZ~bPm  ~h=~bPn 


(4) 


The  pullout  length  ds  is  equal  to  2w,  where  w  is  the  deflection  of 
the  top  or  bottom  beam.  Hence,  the  governing  equation  takes  the 
form: 


El 


d4w 

~dh? 


w 


-2  pm-r=-bpn 


(5) 


The  origin  of  the  x-coordinate  is  assumed  to  be  the  point  where 
the  bridging  zone  begins.  One  should  note  that  the  origin  moves  as 
the  crack  propagates.  The  four  boundary  conditions  (BCs)  for  bridg¬ 
ing  zone  are: 


w(0)  =  w0,  (0  <  w0  <  h/2) 

d2w  d3w 

Mx(0)=apVz(0)  =>  -T(0)  =  op-T(0) 

dx  dx 

w(c)  =  0 


dw 

dx 


(c)  =  0 


(6) 


In  the  above  equation  w0  is  a  prescribed  deflection  at  the  origin  of 
the  coordinate  system  as  shown  in  Fig.  2c.  This  value  will  be  less 
than  h/2  in  the  beginning  and  will  increase  to  a  maximum  value 
of  h/2  as  the  DCB  is  loaded.  Note  that  the  slip  distance  of  the  z- 
pin  is  equal  to  the  total  opening  of  the  DCB  at  that  location,  i.e. 
2w(x)  =  ds(x).  As  the  crack  propagates,  the  bridging  zone  will  also 
move  with  the  crack,  but  the  crack  opening  displacement  will  re¬ 
main  as  h  at  the  origin  with  2w(0)  =  ds( 0)  =  h.  The  terms  V  and  M, 
respectively,  are  the  bending  moment  and  transverse  shear  force 
on  the  beam  cross  section. 

One  should  note  that  the  bridging  length  c  is  still  an  unknown.  It 
can  be  determined  from  the  fact  that  the  strain  energy  release  rate 
at  the  actual  crack  tip  should  be  equal  to  the  mode  I  fracture  tough¬ 
ness  at  the  instant  of  crack  propagation.  The  energy  release  rate 
can  be  determined  from  the  equation  derived  by  Sankar  and  Sonik 
[18]  for  beam-like  specimens  containing  delamination. 

(M(c))2  El(d2w(c))2 

G—bEr-^\-^r)  (7) 


Thus  the  condition  for  determining  c  is: 
El(d2w(c)Y 

b{-*r)  ~Ck 


(8) 


Before  we  solve  the  above  equations  we  will  non-dimensional- 
ize  the  equations  and  BCs  appropriately.  Normalizing  the  length 
dimensions  by  h  and  forces  by  Eh2,  the  governing  equation  and 
the  BCs  take  the  following  form: 
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d4w  __  _ 

^4- -2 PmW=-pm 


w(0)  =  w0,  (0  <  w0  <  1  /2) 

d2w  _  d3w 

dP  (0)  “  Qp  dx3  (0) 
w(c)  =  0 

dw 

dx 

where 


-(c)  =  0 


X  =  T- ,  dp  =~r- 

h  y  h 


-r-C  ,J,  _  W  -jnrl  n  _12P™ 

c  h  ’  w  h  and  Pm  E 


The  equation  for  determining  c  (Eq.  (8))  takes  the  form 


d2w(c) 

~W~ 


=  G/c 


at  0.5  with  additional  increment  of  ap.  This  is  because  the  bridging 
(9)  zone  is  fully  developed  and  it  moves  with  the  crack  tip  as  it 
advances. 

The  procedures  to  solve  the  above  set  of  equations  are  shown  in 
the  flow  chart  depicted  in  Fig.  4.  The  initial  data  includes  the  beam 
properties,  characteristics  of  the  z-pins  and  the  fracture  toughness 
(10)  G/c.  The  deflection  at  the  beginning  of  the  bridging  zone  (x  =  0)  be¬ 
gins  to  increase  as  the  load  is  applied. 

When  w(0)  =  0.5,  bridging  zone  is  fully  developed  and  hereafter 
the  deflection  of  the  beginning  of  bridging  zone  is  constant.  How¬ 
ever,  the  apparent  crack  length  increases  during  crack  propagation. 
We  need  to  use  an  iterative  procedure  as  the  bridging  length  (c)  is 
not  known  a  priori.  The  strain  energy  release  rate  condition  at  the 
right  end  of  the  bridging  zone  (x  =  c)  as  given  by  Eq.  (11)  is  then 
used  to  check  for  correct  value  of  c.  When  a  given  c  satisfies  Eq. 
(11),  then  the  procedure  to  determine  c  is  terminated.  After  bridg¬ 
ing  zone  is  fully  developed,  bridging  length  corresponding  to  every 
^  '  increment  of  crack  length  can  be  determined. 


where  the  non-dimensional  fracture  toughness  is  given  by 

r  -  UGic 
G/c  —  ~Efr 

The  solution  for  the  governing  Eq.  (9)  is: 
w(x)  =  Ci  cos  lx  +  C2  sin  lx  +  C3  cosh  lx  +  C4  sinhlx  +  ^  (12) 

where  1  =  $J2  pm 

The  boundary  condition  at  the  point  x  =  0  varies  since  deflec¬ 
tion  at  this  point  w0  increases  from  zero  at  the  beginning  of  loading 
to  0.5  when  the  bridging  zone  is  completely  developed.  Once  the 
bridging  zone  is  completely  developed  this  value  remains  constant 


Fig.  4.  Flowchart  of  the  procedures  for  solving  the  non-dimensional  governing 
equation. 


3.  Verification 


3  A.  Sample  problem  for  verifying  the  non-dimensional  analytical 
model 


In  order  to  verify  the  analytical  model  a  z-pinned  composites 
designed  by  Cartie  [7]  as  shown  in  Fig.  5  was  selected.  This  z-pin¬ 
ned  composite  was  also  used  in  many  previous  works  [7,15].  The 
configurations  and  material  properties  are  represented  in  Table  2. 
These  values  were  used  for  both  analytical  model  and  FE  simula¬ 
tion  presented  in  the  next  section. 

First,  load-deflection  curve  (Fig.  6)  was  found  for  ap  values 
ranging  from  52.25  mm  to  62.25  mm.  Meanwhile  the  crack  length 
as  the  summation  of  the  apparent  crack  length  (ap)  and  the  bridg¬ 
ing  length  (c)  obtained  from  the  procedure  depicted  in  Fig.  4  was 
computed.  The  variation  of  the  bridging  length  as  a  function  of 
DCB  deflection  3  is  shown  in  Fig.  7.  The  bridging  length  initially  in¬ 
creases  with  loading  until  the  bridging  zone  is  fully  developed, 
where  the  bridging  length  has  the  maximum  value.  In  the  begin¬ 
ning,  the  apparent  crack  length  remains  constant  and  the  crack 
propagation  is  only  due  to  evolution  of  the  bridging  zone.  Once 
the  bridging  zone  is  fully  developed  there  is  little  change  in  the 
bridging  length  thus  the  crack  propagation  is  almost  due  to  the 
apparent  crack  length.  In  fact  there  is  a  slight  decrease  in  the  bridg¬ 
ing  length  after  it  reaches  the  maximum  as  the  bending  moment 
due  to  the  transverse  loads  applied  is  a  function  of  the  apparent 
crack  length.  If  a  pair  of  couple  is  applied  instead  of  a  pair  of  forces 
as  in  the  standard  DCB  test,  then  the  bending  moment  will  remain 
constant  as  the  crack  propagates,  and  one  can  see  a  steady  state 
bridging  length. 

The  force  F  acting  on  the  DCB  and  opening  displacement  <5  at  the 
end  of  the  beam  can  be  obtained  using  the  following  relations: 


F  = -El 


d3w 

d7 


(13) 


3  =  2 


|w(0)  +  ap 


dw 

dx 


x=0 


1  d3w 
3  d7 


x=0 


(14) 


The  variations  of  transverse  force  and  crack  length  with  increments 
in  opening  displacement  are  shown  in  Fig.  6.  Initially  the  transverse 
force  increases  with  opening  displacement  and  begins  to  drop  as 
interlaminar  crack  in  unreinforced  region  propagates.  However, 
the  transverse  force  increases  during  the  development  of  the  bridg¬ 
ing  zone  and  decreases  again  with  the  movement  of  the  fully  devel¬ 
oped  bridging  zone  and  the  new  crack  surface. 
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Fig.  5.  DCB  specimen  reinforced  by  z-pins  [7]. 


Table  2 

The  various  dimensions  and  properties  of  the  DCB  used 
in  the  numerical  simulation  [6]. 


b 

20  mm 

h 

1.6  mm 

E i 

138  Gpa 

e2 

11  GPa 

Vi2 

0.34 

Gl2 

4.4  GPa 

G/c 

258  N/m 

do 

49  mm 

Fm 

18.43  N 

z-Pin  density 

0.5% 

z-Pin  diameter 

0.28  mm 

In  order  to  estimate  the  apparent  fracture  toughness  during 
crack  propagation  we  use  bending  moment  at  the  apparent  crack 
tip  as  a  measure  of  increased  fracture  toughness.  This  is  similar 
to  calculating  the  fracture  toughness  of  unreinforced  beam  using 
Eq.  (7).  Then  the  apparent  fracture  toughness  is  defined  as: 

_  EI  /d2w(0)\ 2 

G,c-app  ~  J  J  (15) 

The  above  relation  can  be  non-dimensionalized  as: 

<w(t)  (16) 


Fig.  7.  Bridging  length  and  crack  length  as  a  function  of  DCB  deflection. 


The  variation  of  apparent  fracture  toughness  during  crack  prop¬ 
agation  is  shown  in  Fig.  8.  The  apparent  fracture  toughness  in¬ 
creases  during  development  of  the  bridging  length  since  the 
bending  moment  required  at  the  beginning  of  bridging  zone  to 


Fig.  6.  Load  and  crack  length  variation  as  a  function  of  DCB  deflection. 
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Fig.  9.  Traction-separation  law  for  the  cohesive  element. 
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Fig.  10.  Cohesive  and  spring  elements  in  the  FE  model  of  the  DCB. 


Fig.  8.  Variation  of  apparent  fracture  toughness  during  crack  propagation. 


overcome  the  bridging  force  is  increases  gradually.  Maximum 
attainable  apparent  fracture  toughness  can  be  predicted  by  the 
relationship  expressed  in  terms  of  interlaminar  fracture  toughness 
and  bridging  force  due  to  z-pins.  The  relationship  is  exactly  the 
same  as  the  equation  based  on  energy  balance.  As  the  crack  prop¬ 
agates  it  has  to  overcome  the  frictional  forces  in  the  z-pins.  The 
amount  of  extra  work  done  is  equal  to  the  area  under  the  load- 
deflection  diagram  in  Fig.  3.  Thus  the  maximum  apparent  fracture 
toughness  can  be  derived  as: 

GiC-app-max  —  G/c  +  hi  0 

Multiplying  throughout  by  (12/Eh)  we  obtain  the  above  relation  in  a 
non-dimensional  form  as: 

GlC-app-max  =  Gfc  ^Pm  (18) 

This  maximum  value  can  be  realized  when  a  pair  of  couple  is 
applied  instead  of  a  pair  of  transverse  forces.  As  shown  in  Fig.  8 
apparent  fracture  toughness  for  transverse  loading  is  always  less 
than  the  maximum  value.  Furthermore,  our  definition  of  apparent 
fracture  toughness  in  Eq.  (16)  is  also  slightly  different  from  the 
energy  based  relation  given  in  Eq.  (18)  making  the  apparent  frac¬ 
ture  toughness  under  transverse  loading  less  than  the  maximum 
realizable  value.  This  is  useful  in  evaluating  the  increase  in  appar¬ 
ent  fracture  toughness  for  a  given  design  of  z-pins.  Moreover  this 
value  is  related  to  strain  exerted  in  the  ligaments  of  DCB  and  will 
be  discussed  later. 

3.2.  FE  simulation  for  verifying  the  non-dimensional  analytical  model 

For  the  sake  of  comparison  with  the  analytical  model,  FE  simu¬ 
lation  of  the  same  specimen  using  the  finite  element  software, 
ABAQUS®,  was  performed.  The  FE  model  containing  discrete  z-pins 
can  verify  the  analytical  model  where  z-pins  are  smeared  and  rep¬ 
resented  by  distributed  traction.  The  three-dimensional  FE  analysis 
was  used  to  simulate  crack  propagation  in  the  DCB  specimen  with 
z-pins.  The  specimen  was  modeled  using  shell  elements  (S4),  and 
cohesive  elements  (COH3D8)  were  used  to  simulate  progressive 
delamination. 

The  behavior  of  cohesive  elements  can  be  characterized  by  bi¬ 
linear  traction  separation  law  (Fig.  9)  [19]  given  by: 


cr  =  (1  -  D)Kd 
where 

r°,  d<d0  n9) 

D=|$8j.  d0<d<df 

1 1,  df  <  d 

where  a  is  traction,  I<  is  stiffness,  D  is  damage  variable,  d  is  displace¬ 
ment,  d0  is  displacement  at  damage  initiation  and  df  is  final  dis¬ 
placement  . 

The  above  parameters  were  taken  as:  I<=  106N/mm, 
am  =  35  Mpa  and  GIC=  0.258  N/mm  [7,19]. 

Nonlinear  spring  element  (CONN3D2)  whose  behavior  is 
defined  by  linearly  decreasing  force  (Fig.  3)  for  each  z-pin  was  also 
implemented  between  the  two  ligaments  of  the  DCB  (Fig.  10). 

The  load-deflection  curves  and  the  variations  of  crack  lengths 
with  deflection  are  shown  in  Fig.  6.  Both  FE  and  analytical  model  re¬ 
sults  are  presented  with  that  from  the  experiment  by  Cartie  [6].  The 
agreement  between  the  analytical  model  and  the  FE  simulations  is 
satisfactory  for  both  load-deflection  and  delamination  length.  The 
slight  discrepancy  in  the  initial  slope  of  the  load-deflection  curve 
between  the  analytical  model  and  the  FE  model  is  due  to  the  fact 
that  the  analytical  model  uses  Euler-Bernoulli  beam  theory 
whereas  FE  model  considered  shear  deformation  as  it  occurs  in  test 
specimens.  The  good  comparison  between  the  analytical  model  and 
the  results  from  FEA  and  experiments  suggests  that  the  discrete 
bridging  force  can  be  represented  by  a  distributed  traction. 

4.  Discussion  of  results  from  the  non-dimensional  model 

4.1.  Parametric  studies  using  the  non-dimensional  analytical  model 

Our  goal  was  to  study  the  effects  of  inherent  interlaminar  frac¬ 
ture  toughness  of  the  composite  material  G/c  and  the  non-dimen¬ 
sional  frictional  force  pm  on  (i)  the  maximum  apparent  fracture 
toughness  G/C_app_max  and  (ii)  steady  state  bridging  length  c.  Such 
relationships  are  extremely  useful  in  design  process  to  evaluate 
the  influence  of  design  variables  on  performance  of  the  composite. 
In  this  parametric  study  we  assume  that  the  DCB  is  loaded  by  end 
couples  instead  of  transverse  forces.  This  assures  steady  state  crack 
propagation  in  the  beam  and  effect  of  increasing^crack  length  on 
the  results  is  thus  eliminated.  In  the  simulations,  Glc  is  varied  from 
10-7  tolO-4  and  pm  ranged  from  10-10  (almost  zero  friction  repre¬ 
senting  the  case  of  unreinforced  beam)  to  0.1. 

Fig.  11  shows^that  the  steady  state  bridging  length  decreases 
with  increasing  Glc  and  pm.  If  the  composite  material  is  inherently 
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Fig.  11.  Non-dimensional  steady  state  bridging  length. 
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Fig.  12.  Apparent  fracture  toughness  as  a  function  of  maximum  friction  force  and  interlaminar  fracture  toughness. 


tough,  it  will  not  allow  a  larger  bridging  length.  Similarly  a  large 
friction  force  also  will  reduce  the  bridging  length.  The  maximum 
apparent  fracture  toughness  (G/C_app_max)  is  another  measure  to 
evaluate  the  effect  of  translaminar  reinforcement.  Fig.  12  shows 
that  the  maximum  apparent  fracture  toughness,  computed  by  Eq. 
(16),  varies  linearly  with  increasing  GIC_app_max  and  pm.  This  result 
is  in  accordance  with^Eq.  (18).  Since  the  range  of  pm  is  wider  than 
that  of  G/c  in  reality,  G/C_app_max  is  more  sensitive  to  the  value  of  pm. 

In  other  words  the  apparent  fracture  toughness  is  dominated  by 
the  translaminar  reinforcement.  From  the  results  shown  in  Figs.  1 1 
and  12,  it  is  clear  that  higher  friction  force  exerted  by  the  z-pins 
provides  enhanced  fracture  toughness  and  at  the  same  time  re¬ 
duces  the  bridging  length.  The  latter  is  important  to  maintain  the 
stiffness  of  the  structure  for  a  larger  bridging  length  leads  to  reduc¬ 
tion  in  the  stiffness  of  the  structure. 


4.2.  Maximum  allowable  translaminar  reinforcement 

Even  though  large  frictional  force  between  the  z-pin  and  the 
surrounding  matrix  material  is  desirable  for  increased  fracture 
toughness,  a  frictional  force  beyond  a  critical  value  will  cause  the 


beam  to  fail.  The  maximum  normal  strain  in  a  beam  cross  section 
is  given  by 

h 


£max  — 


d2w 


dxz 


(20) 


Note  the  strain  is  already  non-dimensional  and  the  right  hand 
side  of  the  above  equation  can  be  written  as: 

1  cTw  K 
2 


dx2 


(21) 


where  k  is  the  non-dimensional  curvature  which  has  the  maximum 
value  at  x  =  0.  Using  Eqs.  (16)  and  (18)  we  obtain 


k=  \/G,c  +  lp„ 


(22) 


Let  us  assume  the  allowable  strain  in  the  composite  is  given  by 


Gu. 


Then, 

G/c  +  2^m  <  u 


(23) 
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(b)  Force  Input 


(b)  Moment  Input 

Fig.  13.  Comparison  of  strain  with  different  loading  types. 


The  above  equation  provides  an  upper  limit  on  the  z-pin  density 
which  should  be  taken  into  consideration  in  the  design  of  transla- 
minar  reinforcements.  Thus  the  allowable  density  or  maximum 
frictional  force  of  z-pins  can  be  determined  at  a  given  geometry 
and  material  properties  of  a  composite  material. 

In  order  to  verify  the  result  for  the  maximum  allowable  transla- 
minar  reinforcement  given  in  Eq.  (24),  we  performed  finite  element 
analysis  of  a  specimen.  The  properties  used  were:  Glc  =  0.258 
N/mm,h  =  1.6  mm,E  =  138  GPa,N/m  =  1.58  N/mm2  [7],  The  FE 
simulation  was  used  to  calculate  the  maximum  strain  in  the  liga¬ 
ments  of  the  DCB  specimen  at  the  instant  of  crack  propagation, 
and  it  was  compared  with  the  eu  obtained  from  Eq.  (24).  Two 
loading  cases,  namely  end  couples  and  transverse  forces,  were 
considered.  For  the  set  of  properties  used  Eq.  (24)  yields  a  maxi¬ 
mum  strain  eu  =  4.5  x  10-3.  Note  that  if  Nfm  >  1.58  N/mm2  then 
the  maximum  strain  exceeds  the  ultimate  strain.  From  the  FE 
simulations  we  obtained  £max  =  4.521  x  10-3  for  the  case  of  end 
couples  and  emax  =  3.906  x  10-3  for  the  case  of  transverse  loading 
of  the  DCB  (Fig.  13).  The  results  indicate  the  strains  in  the  beam 
are  higher  for  the  moment  loading  compared  to  that  due  to  the 
transverse  loading.  However,  in  a  practical  structure  the  delamina¬ 
tion  will  experience  a  combination  of  shear  and  moment  loadings, 
and  hence  the  conservative  value  should  be  used.  That  means  the 
maximum  allowable  stitch  density  given  in  Eq.  (24)  should  be  used 
although  it  is  applicable  only  to  moment  loading  case,  and  hence 
conservative. 


E,  Ell 

SNEG,  (fraction  = 
(Avg:  75%) 

+4.521e-03 
+3. 765e-03 
+3.008e-03 
+2. 252e-03 
+1. 495e-03 
+7.390e-04 
-1. 743e-05 
-7.739e-04 
-1.530e-03 
-2.287e-03 
-3.043e-03 
-3.800e-03 
-4.556e-03 


-1.0) 


Fig.  14.  Tri-linear  bridging  law  is  indicated  by  solid  line.  The  dotted  line  is  the 
linear  softening  law  in  Fig.  3.  The  areas  under  the  force-displacement  diagram  for 
both  laws  are  the  same. 

From  the  above  equation  one  can  derive 

Pm  <  2(4 £2  —  G/c) 

or  (24) 

Pm  =  Nfm<lEs2u 


4.3.  Effect  of  tri-linear  bridging  law 

So  far  we  have  used  the  linear  softening  bridging  law  for  devel¬ 
oping  the  non-dimensional  analytical  model  and  verifying  the 
same  using  FEA.  Although  this  bridging  law  is  simple,  it  may  not 
be  realistic.  Dai  et  al.  [6]  performed  pull-out  tests  using  z-pins  of 
various  diameters  to  determine  the  actual  bridging  law.  They 
found  that  a  high  value  of  debonding  force  was  reached  before 
the  debonding  of  the  pins  began.  After  the  debonding  was  initiated, 
the  pull-out  force  dropped  to  a  lower  value  before  reducing  line¬ 
arly  to  zero  value  as  the  pins  were  pulled  out  steadily  against  fric¬ 
tional  force  that  seemed  to  have  a  constant  coefficient  of  friction. 
They  represented  this  pullout  behavior  by  a  tri-linear  bridging 
law  as  depicted  in  Fig.  14.  The  tri-linear  behavior  can  be  attributed 
to  an  elastic  deformation  (from  (0,  0)  to  (0.0185,35.3)  in  Fig.  14), 
debonding  from  surrounding  matrix  (from  (0.0185,35.3)  to 
(0.17, 14.86)  in  Fig.  14)  and  slip-out  of  z-pins  (from  (0.17, 14.86) 
to  (1.6,0)  in  Fig.  14). 

The  finite  element  simulations  of  the  DCB  described  in  Section 
3.2  were  repeated  with  the  tri-linear  bridging  law.  We  did  not  at- 
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Fig.  15.  Load-deflection  curve  from  various  bridging  laws. 
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tempt  analytical  solution  as  it  will  be  very  complicated  due  to  the 
piecewise  continuous  bridging  law.  The  numerical  values  for  the 
tri-linear  law  were  chosen  such  that  the  area  under  the  force-dis¬ 
placement  curve,  which  is  related  to  the  increase  in  fracture  tough¬ 
ness,  is  equal  to  the  linear  softening  law  used  in  the  earlier 
example.  That  is  the  areas  under  the  two  force-displacement 
curves  in  Fig.  14  are  equal  to  each  other.  The  resulting  load-deflec¬ 
tion  diagram  of  the  z-pinned  DCB  specimen  is  shown  in  Fig.  15.  The 
load  deflection  behavior  is  almost  identical  to  that  obtained  using 
the  linear  bridging  law  indicating  that  the  apparent  fracture  tough¬ 
ness  is  same  in  both  cases.  It  must  be  noted  that  the  peak  force  in 
the  tri-linear  case  is  almost  twice  as  that  of  the  linear  softening 
law.  The  variations  of  crack  length  as  a  function  of  the  DCB  speci¬ 
men  opening  (deflection)  are  also  similar  in  both  cases.  This  indi¬ 
cates  that  the  details  of  the  bridging  law  do  not  seem  to  affect 
the  global  behavior  of  the  specimen  as  long  as  the  energy  dissi¬ 
pated  by  the  pins  is  properly  accounted  for. 

5.  Summary  and  conclusions 

Mode  I  delamination  propagation  in  DCB  specimens  containing 
z-pins  is  studied.  A  simple  analytical  model  based  on  linear 
softening  type  bridging  law  for  the  z-pins  has  been  developed  and 
suitable  non-dimensional  parameters  have  been  identified.  The 
load-deflection  curve  of  the  DCB  specimen  was  calculated  using 
the  analytical  model.  It  is  seen  that  the  bridging  zone,  wherein 
the  pins  are  partially  pulled  out,  develops  as  the  crack  propagates, 
but  attains  a  steady  state  value.  The  length  of  the  bridging  zone  is 
a  function  of  the  mode  I  fracture  toughness  and  the  frictional  force 
between  the  z-pins  and  the  surrounding  material.  An  expression 
was  derived  for  the  apparent  or  effective  fracture  toughness  values. 
Although  increase  in  frictional  force  as  the  z-pins  increases  the 
fracture  toughness,  there  is  an  upper  limit  to  this  friction  as  the 
DCB  ligaments  would  break  if  the  friction  is  very  high.  The  limiting 
value  of  the  pin  friction  is  derived. 

The  efficacy  of  the  analytical  model  was  evaluated  by  the  sim¬ 
ulation  of  the  DCB  specimen  using  finite  element  simulations.  In 
the  FE  model  the  delamination  propagation  was  simulated  by 
cohesive  elements  and  the  z-pins  were  modeled  as  discrete  nonlin¬ 
ear  elements.  The  results  for  load-deflection  curve  and  the  crack 
bridging  zone  length  agreed  quite  well  with  the  analytical  model. 
As  an  alternative  to  the  linear  bridging  law,  a  more  realistic  tri-lin¬ 
ear  bridging  law  was  used  in  the  FE  simulations.  It  is  found  that  the 
global  delamination  behavior  of  the  specimen  was  not  affected 
much  as  long  as  the  energy  dissipated  by  the  pins  is  kept  the  same. 

The  non-dimensional  model  with  few  parameters  will  serve  as  a 
design  tool  when  translaminar  reinforcements  such  as  z-pins  are 
selected  for  laminated  composite  structures  in  order  to  improve 
their  fracture  toughness.  The  analytical  models  will  also  be  useful 


in  optimization  studies  and  simulation  of  large  composite  struc¬ 
tures  containing  translaminar  reinforcements. 
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